#' ---
#' title: "Descriptive Statistics"
#' author: "Lukas F. Stoetzer"
#' date: "Feb 2022"
#' ---

  
# Read Data
  datw <- readRDS("dat_wide.RDS")
  
# Party Lab
  party_lab <- read_csv("dat_partylab.csv")
  
# Panel characteristics + Quota for each Country ============
  
  # Prepare Target Quota
  quota <- read_csv("international_quotas.csv") 
  
  # Prepare Data Frame        
  quata_tab <- datw %>% 
    mutate(Male = ifelse(gender == 1, 1,0),
           Female = ifelse(gender == 1, 1,0), 
           age = 2020 - as.numeric(as.character(as_factor(datw$birthyear))),
           Age1829 = ifelse(age %in% 18:29, 1,0),
           Age3039 = ifelse(age %in% 30:39, 1,0),
           Age4049 = ifelse(age %in% 40:49, 1,0),
           Age5059 = ifelse(age %in% 50:59, 1,0),
           Age6075 = ifelse(age %in% 60:75, 1,0), 
           edulow = case_when(education %in% 1:2 & cntry_resp %in% c("ger","usa","bra","pol") ~ 1,
                              education %in% 1 & cntry_resp %in% c("ita") ~ 1,
                              TRUE ~ 0),
           edumid = case_when(education %in% 3 & cntry_resp %in% c("ger","usa") ~ 1,
                              education %in% 2 & cntry_resp %in% c("ita") ~ 1,
                              education %in% 3:4 & cntry_resp %in% c("bra","pol") ~ 1,
                              TRUE ~ 0),
           eduhigh = case_when(education %in% 4:6 & cntry_resp %in% c("ger","usa") ~ 1,
                               education %in% 4:6 & cntry_resp %in% c("ita") ~ 1,
                               education %in% 5:6 & cntry_resp%in% c("bra","pol") ~ 1,
                               TRUE ~ 0)) %>% 
    select(cntry_resp,Male, Female, starts_with("Age"), edulow, edumid, eduhigh,-age)  %>%
    group_by(cntry_resp) %>%
    summarise_all(mean) %>%
    gather(var, val, -cntry_resp) %>% 
    spread(cntry_resp, val) %>% 
    gather(cntry,Share_Survey,-var) %>%
    mutate(cntry = case_when(
      cntry == "bra" ~ "Brazil",
      cntry == "ger" ~ "Germany",
      cntry == "ita" ~ "Italy",
      cntry == "pol" ~ "Poland",
      cntry == "usa" ~ "USA"
    )) %>%
    left_join(.,gather(quota,cntry,Quota,-var)) %>%   
    mutate(var = case_when(
      var == "Age1829" ~ "Age 18-29",
      var == "Age3039" ~ "Age 30-39",
      var == "Age4049" ~ "Age 40-49",
      var == "Age5059" ~ "Age 50-59",
      var == "Age6075" ~ "Age 60-75",
      var == "eduhigh"  ~ "Education High",
      var == "edulow"   ~ "Education Low",
      var == "edumid"   ~ "Education Middle",
      TRUE ~ var
    ))
    

  # Make Table
  select(quata_tab,"Variable" = "var", "Share Survey" = "Share_Survey","Quota") %>% 
  kable(.,format = "latex", digits = 2, booktabs=T,longtable = T,
        caption = "Quota and survey shares in different countries ") %>%
    kable_styling(latex_options = "hold_position") %>% 
    pack_rows("Brazil", 1, 10, latex_gap_space = "1em") %>%
    pack_rows("Germany", 11, 20, latex_gap_space = "1em") %>%
    pack_rows("Italy", 21, 30, latex_gap_space = "1em") %>%
    pack_rows("Poland", 31, 40, latex_gap_space = "1em") %>%
    pack_rows("USA", 41, 50, latex_gap_space = "1em")   %>% 
    cat(.,file="appendix_table4.tex")

# Attention check share country-wise ===============
  
  datw %>% 
    group_by(cntry_resp) %>% 
    summarise("Passed Attention Check" = mean(cjattcheck,na.rm=T)) %>% 
    mutate(cntry_resp = case_when(cntry_resp == "bra" ~ "Brazil",
                                  cntry_resp == "ger" ~ "Germany",
                                  cntry_resp == "ita" ~ "Italy",
                                  cntry_resp == "pol" ~ "Poland",
                                  cntry_resp == "usa" ~ "USA")) %>%
      rename("Country" = cntry_resp)   %>%
      kable(format = "latex", digits = 2, booktabs=T) %>%
      kable_styling(latex_options = "striped") %>% 
      cat(.,file="appendix_table5.tex")
    

  
# Affective Polarization between Party Options ===================
  
  # Reshape Data
  dat_plot <- datw %>%
    select("cntry"=cntry_resp, psid,starts_with("partyfeeling")) %>%
    mutate(affec_pol = case_when(
      cntry == "usa" ~ partyfeeling_1 - partyfeeling_2,
      cntry == "ger" ~ partyfeeling_7 - partyfeeling_4,
      cntry == "ita" ~ partyfeeling_2 - partyfeeling_5,
      cntry == "pol" ~ partyfeeling_1 - partyfeeling_2,
      cntry == "bra" ~ partyfeeling_2 - partyfeeling_1)) %>%
    mutate(cntry = case_when(cntry == "bra" ~ "Brazil",
                             cntry == "ita" ~ "Italy",
                             cntry == "ger" ~ "Germany",
                             cntry == "usa" ~ "USA",
                             cntry == "pol" ~ "Poland"))


  ggplot(dat_plot) +
    geom_density_ridges_gradient(aes(x = affec_pol, y = cntry)) +
    scale_fill_manual(
      name = "Affective Polarization", values = c("#FF0000A0", "#A0A0A0A0", "#0000FFA0"),
      labels = c("Left", "Middle", "Right")
    ) +
    geom_vline(xintercept=c(-5,5), alpha=0.3) +
    xlab("Affective Polarization (Difference in Feeling Thermometer Scores)") +
    ylab("") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()


  # Save
  ggsave("appendix_figure7.pdf", width = 12, height=9)
  
  
  ggplot(dat_plot) +
    geom_density_ridges_gradient(aes(x = affec_pol, y = cntry,fill=(..quantile..)),
                                 quantile_lines = c(TRUE), 
                                 quantiles =c(0.33,0.66)) +
    scale_fill_manual(
      name = "Affective Polarization", values = c("#FF0000A0", "#A0A0A0A0", "#0000FFA0"),
      labels = c("Left", "Middle", "Right")
    ) +
    xlab("Affective Polarization (Difference in Feeling Thermometer Scores)") +
    ylab("") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()
  
  # Save
  ggsave("appendix_figure14.pdf", width = 12, height=9)
  
  

# Job Relevance in different countries ===================
  
  # Jab Labels
  item_lab <- data.frame("job_item"=1:6,
             "job" =c("Doctors","Nurses","Cooks","Journalists",
                      "Professors","Police Officers") )


  # Prepare Data
  d_job <- datw %>% 
    select(cntry_resp, starts_with("jobs_relevance")) %>%
    gather(job, val, -cntry_resp) %>%
    separate(job,c(NA,NA,"job_item")) %>%
    mutate(job_item = as.numeric(job_item)) %>%
    bind_rows(.,mutate(.,cntry_resp = "pooled")) %>%
    group_by(cntry_resp,job_item) %>%
    summarise("Percent_Very_Essential" = mean(val==3,na.rm=T)) %>%
    left_join(.,item_lab) %>%
    arrange(cntry_resp,-Percent_Very_Essential) %>%
    ungroup() 

  
  # Figure 
  d_job <- d_job %>% 
    mutate(job = factor(job, levels = rev(c("Doctors","Nurses", "Police Officers","Professors","Cooks","Journalists"))),
           cntry_resp = factor(cntry_resp, 
                               levels = c("pooled","usa", "ger","ita","pol","bra"),
                               labels = c("Countries combined","USA","Germany","Italy","Poland","Brazil")))
  
  ggplot(d_job,aes(y=job,x=Percent_Very_Essential)) +
    geom_point() +
    facet_wrap(~ cntry_resp) +
    xlab("Very essential job (%)") +
    ylab("") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()
 
  ggsave("appendix_figure6.pdf",width = 9, height = 6)
  
# Experience with COVID ==============
  
  # Aggreagte Data
  dat_expCOVID  <- datw %>%
    select(cntry_resp, corona_infected) %>%
    mutate(corona_infected = as_factor(corona_infected)) %>%
    bind_rows(.,mutate(.,cntry_resp = "pooled")) %>%
    na.omit() %>%
    group_by(cntry_resp,corona_infected) %>%
    summarise(n = n()) %>%
    mutate(perct = n / sum(n), 
           cntry_resp = factor(cntry_resp, 
                               levels = c("pooled","usa", "ger","ita","pol","bra"),
                               labels = c("Countries combined","USA","Germany","Italy","Poland","Brazil")))
  
  # Plot 
  ggplot(dat_expCOVID) +
    geom_bar(aes(x=corona_infected,y=perct),stat="identity") +
    facet_wrap(~ cntry_resp) +
    xlab("Know People infected with COVID") +
    ylab("Percentage") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()
  
  ggsave("appendix_figure8.pdf",width = 12, height = 9)
  
# Risk Battery & Compliance by partisan groups ===============
  
  # Prepare Data
  data_risk <- datw %>% 
    bind_rows(., mutate(.,cntry_resp="pooled")) %>%
    group_by(cntry_resp,resp_pid) %>%
    summarise("mean" = mean(riskybehav_score,na.rm=T),
              "se" = sd(riskybehav_score,na.rm=T)/ sqrt(n()-1)) %>%
    mutate(cntry_resp = factor(cntry_resp, 
                               levels = c("pooled","usa", "ger","ita","pol","bra"),
                               labels = c("Countries combined","USA","Germany","Italy","Poland","Brazil")))
  

  # Plot
  ggplot(data_risk) +
    geom_pointrange(aes(x=resp_pid, 
                        y= mean,
                        ymin= mean - 2*se,
                        ymax = mean + 2*se)) +
    coord_flip() +
    facet_wrap(~ cntry_resp) +
    xlab("Risk behavior index") +
    ylab("Partisan respondent") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()
  
  ggsave("appendix_figure10.pdf",width = 9, height = 6)
  
  
  # Compliance
  dat_compl <- datw %>%
    mutate(corona_precautions_always_follow = ifelse(corona_precautions == 1,1,0)) %>%
    bind_rows(., mutate(.,cntry_resp="pooled")) %>%
    group_by(cntry_resp,resp_pid) %>%
    summarise("mean" = mean(corona_precautions_always_follow,na.rm=T),
              "se" = sd(corona_precautions_always_follow,na.rm=T)/ sqrt(n()-1)) %>%
    mutate(cntry_resp = factor(cntry_resp, 
                               levels = c("pooled","usa", "ger","ita","pol","bra"),
                               labels = c("Countries combined","USA","Germany","Italy","Poland","Brazil")))
  
  
  ggplot(dat_compl) +
    geom_pointrange(aes(x=resp_pid, 
                        y= mean,
                        ymin= mean - 2*se,
                        ymax = mean + 2*se)) +
    coord_flip() +
    facet_wrap(~ cntry_resp) +
    ylab("Percentage always comply with recommendations") +
    xlab("Partisan Respondent") +
    theme(text = element_text(size=12),
          panel.grid.major.y =  element_blank(),
          panel.grid.minor.y =  element_blank()) +
    theme_bw()
    
  ggsave("appendix_figure9.pdf",width = 9, height = 6)
  

    